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An electromechanical Ising machine 
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Solving intractable mathematical problems in simulators composed of atoms, ions, photons or elec¬ 
trons has recently emerged as a subject of intense interest. Here we extend this concept to phonons 
that are localised in spectrally pure resonances in an electromechanical system which enables their 
interactions to be exquisitely fashioned via electrical means. We harness this platform to emulate 
the Ising Hamiltonian whose spin 1/2 particles are replicated by the phase bistable vibrations from 
a parametric resonance where multiple resonances play the role of a spin bath. The coupling be¬ 
tween the mechanical pseudo spins is created by generating thermomechanical two-mode squeezed 
states which impart correlations between resonances that can imitate a ferromagnetic, random or 
an anti-ferromagnetic state on demand. These results suggest an electromechanical simulator for 
the Ising Hamiltonian could be built with a large number of spins with multiple degrees of coupling, 
a task that would overwhelm a conventional computer. 


Electro/optomechanical systems consisting of a high 
quality mechanical resonance embedded in an electri¬ 
cal/optical transduction circuit [THil have emerged as 
versatile platform in which sensors with unprecedented 
sensitivities can be developed m\, in which mechan¬ 
ical nonlinearities can be dynamically engineered and 
harvested HHin] and even quantum mechanics within a 
macroscopic context can be studied HIHIS]- However 
this unprecedented success has had the unfortunate con¬ 
sequence of stifling further advancement of mechanical 
systems beyond these extensively explored concepts. 

Meanwhile quantum simulators have been proposed 
as devices which could be used to solve problems in 
quantum physics that are beyond the capacity of clas¬ 
sical computers El- This notion can also be applied to 
mathematical problems for example the Ising Hamilto¬ 
nian used to describe a spin glass which is characterised 
as an NP-hard problem that cannot be solved via conven¬ 
tional computing protocols m- However in spite of this 
formidable challenge, the ground state spin configuration 
of the Ising Hamiltonian could yield solutions to optimi¬ 
sation problems that can be mapped on to its spin lattice 
and thus its efficient extraction is highly desired [16] . One 
approach to rapidly determining the ground state of the 
Ising Hamiltonian is to employ quantum annealing where 
quantum tunnelling is harnessed to search the energy 
landscape corresponding to a given problem programmed 
into the underlying spin configuration pYHlQ] . Recently 
an alternative and apparently classical approach to this 
problem has emerged where a time multiplexed optical 
parametric resonator network is programmed with an 
Ising Hamiltonian composed of four spins. The ground 
state is then determined by simply activating the net¬ 
work which preferentially resonates in its global poten¬ 
tial minima [20l [21] . Here this concept is developed in a 
frequency multiplexed electromechanical parametric res¬ 
onator and we show that this platform can readily be 
extended to multiple electromechanical parametric res¬ 
onators with arbitrary degrees of coupling namely all the 
necessary prerequisites to solving the Ising Hamiltonian 



FIG. I: The electromechanical system, a, A false colour 
electron micrograph of the coupled mechanical resonators sus¬ 
taining the symmetric and asymmetric vibration modes. The 
measurements were performed at room temperature and in 
a high vacuum and the mechanical vibrations were detected 
via a laser interferometer which was demodulated either in 
a spectrum analyser (SA) or in a phase sensitive detector 
utilising the heterodyne mixing setup detailed in the circuit 
schematic. Five signal generators were employed with two 
piezoelectrically activating the parametric resonances at 2ujn, 
two generating the reference signals for the phase sensitive de¬ 
tectors at (jJn and the fifth creating the coupling between the 
parametric resonances via non-degenerate parametric down- 
conversion when activated at . 


in a regime that has no analytical solution and where its 
numerical solution would demand exorbitant computa¬ 
tional resources. 

The Ising model conceived using the tools of statis¬ 
tical physics describes ferromagnetism and in the ab¬ 
sence of a magnetic field its Hamiltonian is given by 

N 

— Jik<^iO'k with N particles each with two spin states 
ik 
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cr = zbl and the coupling between the and parti¬ 
cles is parametrised by Jik- In the first steps to building 
an electromechanical Ising machine, the fundamentals of 
the Ising Hamiltonian need to be replicated in the elec¬ 
tromechanical domain namely a mechanical pseudo spin 
to encode a, multiple mechanical spins playing the role of 
an N particle bath and finally coupling between the me¬ 
chanical spins Jik whose magnitude and polarity can be 
controlled and has the potential to be extended beyond 
nearest neighbour spins. 

The prototype electromechanical system in which these 
concepts are investigated is shown in Fig. 1 and it con¬ 
sists of two strongly coupled mechanical beams which 
sustain a symmetric (S) and an asymmetric (A) mode 
at Co’g/27r ^245 kHz and ^260 kHz with qual¬ 

ity factors (Q) of 1300 and 2300 respectively (see Fig. 
SI in Supplementary Information (SI)). The mechani¬ 
cal elements are integrated with piezoelectric transduc¬ 
ers which enable the underlying harmonic potential of 
both modes to be parametrically modulated leading to a 
system Hamiltonian given by 

A 

^ (p»/^ + ^n9^/2(l - 2r„ cos(2w„i) 

n=S 

+/3n9^/2)) + AQsQa cos{uj^t + (p) (1) 

where the summation expresses the kinetic and poten¬ 
tial energies from both modes in terms of their position 
Qn and momentum pn [H [22l [23]. The potential en¬ 
ergy term contains three contributions, the second arising 
from the periodic modulation of the mechanical spring 
constant with amplitude F^ at twice the natural mode 
frequency which yields degenerate parametric amplifica¬ 
tion and parametric resonance [24l [25] and the third due 
to the well-known anharmonicity from the Duffing non¬ 
linearity /3n which appears at large displacements [26] . 
A unique feature of the parametric resonance is that it 
can only vibrate with two phases and thus it provides the 
ideal construct within which a classical spin can be har¬ 
boured lanni- The last term describes non-degenerate 
parametric down-conversion when the periodic modula¬ 
tion is activated with a pump {V) amplitude A at the 
sum frequency of both modes ) which re¬ 

sults in the symmetric and asymmetric modes becoming 
amplified and correlated yielding a two-mode squeezed 
state [22] . 

This Hamiltonian can be transformed in the rotating 
frame approximation with the introduction of a dimen¬ 
sionless canonical position coordinate Qn and conjugate 
momentum P^, as detailed in SI, which yields 


a 


E 

Q. 


c 


E 

Q. 


CL 


Qa (pm) 


2 



-1 


-2 - 

- 2-1012 
Qa (pm) 



b 


E 


d 


50 

25 



I 


-50 

- 3 - 2-10123 
Qa (nm) 



) y^^^^^^^^^J 
0.0 0.1 0.2 0.3 0.4 0.5 


Vp(®p) (V,,3) 



y^{2(D^) and V^(2a)^) 


FIG. 2: Mechanical spins and spin coupling, a and b, 

The occupation probabilities of the parametric resonances in 
phase space are measured with the parametric actuation pe¬ 
riodically activated with F^(2cj^) = 0.7 Yrms and V^(2u;^) 
= 0.6 Yrms respectively. The resultant phase portraits in¬ 
dicate that when a parametric resonance is activated, each 
mode evolves from its stationary state at the origin, which is 
broadened by thermomechanical noise fluctuations, to one of 
the two available oscillating states with a tt phase separation 
that can be exploited to host classical spin information, c, 
The phase portrait reconstructed from the cross-quadratures 
of the symmetric and asymmetric mode as function of pump 
intensity, d, The corresponding correlation coefficient reveals 
that the two modes become indistinguishable as the pump am¬ 
plitude is increased where the colours indicate the amplitudes 
in the previous figure, e, A two-mode squeezed state gener¬ 
ated with V^{2ujj,) = 0.2 Yrms which is then enhanced via 
degenerate parametric amplification of both modes, f, The 
corresponding correlation coefficient confirms that the initial 
weakly correlated state can be enhanced by the degenerate 
parametric amplification where the colour coding indicates 
the amplitudes in the previous figure. 


A 


= + JsA<^s<^A 

n=S 


where the first term describes the energy correspond¬ 
ing to the two stable oscillation phases of a paramet¬ 
ric resonance at = 0 and Qn = (Jn\/^^n/^/3n 
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EH W\ and the second term quantifies the coupling 
= A cos(0) uj^ Pg /3^ between the modes 

with their two phases encoding = ±1. Thus this elec¬ 
tromechanical system can be formally mapped onto the 
Ising model composed of two classical spins arising from 
the phase bistable parametric resonances of both modes 
that can be coupled from the non-degenerate parametric 
down-conversion. Pleasingly, this analysis also reveals 
that the magnitude of the coupling between the mechan¬ 
ical spins created by A can be enhanced by the degenerate 
parametric resonances that are activated when T^ > Tn, 
where the damping rate 7^ = I Qn [S], and its polarity 

i.e. the sign of can be selected by simply tuning the 
pump phase 0 [28] . 

To experimentally confirm these expectations, the 
spring constant of either mode is modulated at Vni^uJn) 
via the stress induced by the piezoelectric transduc¬ 
ers namely T^ in the above equations is activated [9]. 
This results in a given modes thermomechanical fluctu¬ 
ations initially being enhanced corresponding to degen¬ 
erate parametric amplification [24] and once the rate at 
which the phonons are generated begins to exceed their 
rate of decay namely T^ > 7 ^ it results in infinite am¬ 
plification corresponding to a parametric resonance as 
depicted in Figs. S 2 |9l[25]. The parametric resonance, 
under periodic driving, is then projected into phase space 
by mixing with a local oscillator locked onto its resonance 
frequency and demodulated in a phase sensitive detec¬ 
tor (see circuit schematic in Fig. 1). This measurement 
yields two time series for the in-phase Qn and quadra¬ 
ture Pn components that enables a phase portrait to be 
reconstructed as shown in Figs. 2a and 2b which prob¬ 
abilistically unveils the parametric resonances from both 
modes being able to vibrate with only two phases sep¬ 
arated by TT radians as described above [9]. This phase 
bistable vibration enables a classical spin to be mimicked 
in the mechanical domain where hence forth the positive 
in-phase component is defined as spin up i.e. = 1 and 
vice versa. 

In order to create coupling between the pseudo spins 
encapsulated by the symmetric and asymmetric modes, a 
mechanical two-mode squeezed state is created by pump¬ 
ing the spring constant of both modes at ) corre¬ 

sponding to A in the above equations [ 22 ] . This results in 
the thermomechanical fluctuations of both modes being 
simultaneously amplified from non-degenerate paramet¬ 
ric down-conversion as shown in Fig. S 2 and projecting 
this output in phase space reveals it to be phase insen¬ 
sitive as shown Fig. S4 [29]. In contrast, the degenerate 
parametric amplification detailed above is phase sensi¬ 
tive with only the in-phase component being amplified 
whilst the quadrature component is squeezed as shown 
in Fig. S3 [24]. However, the simultaneous generation 
of phonons in both modes also results in their vibrations 
becoming correlated which can be identified by recon¬ 
structing their cross-quadratures in phase space namely 


the in-phase component of the asymmetric mode versus 
the quadrature component of the symmetric mode (or 
vice versa) as shown in Fig. 2c (Fig. S4c). The resul¬ 
tant squashed distribution implies that the motion from 
both modes has become intertwined which can statisti¬ 
cally be confirmed by evaluating their correlation coeffi¬ 
cient cov{Q^Pg)/(Q^(p^ where the numerator describes 
the covariance and ( is the standard deviation. The re¬ 
sults of this analysis (also including Qg Pa whose phase 
portraits are shown in Fig. S4) as a function of pump in¬ 
tensity are shown in Fig. 2d which reveals that in this 
process the correlation coefficient tends to -1 indicating a 
perfect anti-correlation that is phonons pairs are simulta¬ 
neously generated in both modes and with an increasing 
rate which renders their vibrations classically indistin¬ 
guishable. 

To transfer this correlation into the parametric res¬ 
onances, a weakly correlated state is first created at 
H^(co’^)=200 mV^^s as shown in Figs. 2e and 2 f. Next 
degenerate parametric amplification is continuously ac¬ 
tivated in both modes and the resultant correlation co¬ 
efficient is monitored as function of both Vg{2ujg) and 
Vj^{2uj^). As expected from equation 2 , the correlations 
between the two modes can be enhanced from the degen¬ 
erate parametric amplification which can readily be seen 
by the narrower squashed distribution in phase space as 
shown in Fig. 2e and quantitatively confirmed via the 
correlation coefficient in Fig. 2 f. Indeed once maximum 
degenerate parametric amplification has been achieved 
and both modes start to parametrically resonate, a per¬ 
fectly correlated state is created even with weak pump 
intensity. Although experimentally the evolution of the 
correlation coefficient across the transition from degen¬ 
erate parametric amplification to parametric resonance 
cannot be seamlessly measured as the phase sensitive 
detectors became overloaded, this transition is readily 
confirmed both numerically and analytically by solving 
equation 1 in this configuration as shown in Figs. S5 and 
S 6 . 

Based on the above results, the essential features of 
the Ising Hamiltonian between two spin 1/2 particles in 
the electromechanical domain can now be demonstrated 
using the pulse sequence depicted in Fig. 3a. First a 
two-mode squeezed state is created which is undetectable 
via thermo mechanical fluctuations namely Vp(uj^) 0 . 
Next, the symmetric mode is activated off-resonance i.e. 
2(co’g + A) which results in no vibration. This is then 
followed by the asymmetric mode being parametrically 
activated at 2uj^ where the tension created by its mo¬ 
tion modifies the spring constant and shifts the natu¬ 
ral frequency of the symmetric mode by A resulting in 
its parametric resonance being simultaneously activated 
as shown in Fig. 3b. In this formation, periodically 
triggering the asymmetric mode results in the symmet¬ 
ric mode becoming activating with the same period and 
monitoring the in-phase component from both modes on 
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FIG. 3: Emulating the Ising model in the mechanical domain, a, The pulse sequence used to demonstrate the funda¬ 
mentals of the Ising Hamiltonian in an electromechanical system, b, The spectral response of the symmetric and asymmetric 
modes in response to the above pulse sequence, c, The temporal response of the mechanical spins to the above pulse sequence 
which are readout via the in-phase component of their motion with + A)) = 0.7 Vrms, ^a(2^a) — Vrms and 

Vj, (ct;p)=600 /jMrms which reveals ferromagnetic ordering with the two spins always in parallel alignment where the spin orien¬ 
tations are visualised by the arrows with red (green) indicating the symmetric (asymmetric) mode. All the temporal outputs are 
measured over 1000 seconds to enable statistically reliable analysis, d, The correlation coefficient extracted from the switching 
outputs as a function of pump intensity. At the largest pump amplitude the mechanical spins exhibit ferromagnetic ordering 
but as the pump is reduced to zero, the spins become disordered undergoing an effective second order phase transition with 
the pump playing the role of temperature and a critical temperature (Tc) of 500 fiVrms- e, A segment of the temporal output 
used to extract the above correlation coefficient with the pump deactivated which clearly reveals that the mechanical spins 
have become decoupled yielding uncorrelated switching events, g, The correlation coefficient extracted from the mechanical 
system in response to the above pulse sequence as a function of the pump phase with (ct;p)=600 fiYrms- As the pump 
phase is adjusted, the mechanical spins exhibiting ferromagnetic ordering at 0 — 0o=67° detailed in the temporal output in f 
transition through a disordered state with random spin alignment to finally an anti-correlated state at 0 — 0o=-63° with the 
spins exhibiting anti-ferromagnetic ordering as detailed in the temporal output in h. Note that the value of 0o is determined 
by the spectral position of the readout probe with respect to the resonances of the two modes. 


resonance enables identification of their spin orientation 
(see Figs. 2a and 2b). Correlations from the two-mode 
squeezing can then be statistically quantified by moni¬ 
toring the resultant switching outputs from both modes. 
For instance, activating the asymmetric mode will result 
in either a spin up or a spin down the choice of which 
is determined by its thermal fluctuations. However, if a 
two-mode squeezed state is created then the switching 
output observed from the symmetric mode will be corre¬ 
lated with that from the asymmetric mode. 

Implementing this protocol with a weak pump yields 
the output shown in Fig. 3c where an up spin from the 
asymmetric mode coincides with an up spin in the sym¬ 


metric mode and vice versa. However if the pump is de¬ 
activated, the two modes become uncoupled and switch 
randomly as shown in Fig. 3e. The correlation coefficient 
can then be extracted from these switching outputs as a 
function of pump intensity as depicted in Fig. 3d. This 
reveals the mechanical pseudo spins can achieve parallel 
alignment when the pump is activated and the corre¬ 
sponding correlation coefficient tends to unity. This par¬ 
allel configuration can be controllably deactivated via the 
pump and the correlation coefficient smoothly transitions 
to zero mimicking an effective second order ferromagnetic 
phase transition. Consequently ) can convincingly 

play the role of Jik in the Ising Hamiltonian. 
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Now in order to control the sign of the coupling, the 
pump phase is adjusted as depicted in Fig. 1 and detailed 
in equation 2. Indeed measuring the two-mode squeezing 
via thermomechanical noise fluctuations in phase-space 
as a function of (j) clearly reveals that the anti-correlated 
state depicted in Fig. 2c can be tuned to a correlated 
state as shown in Fig. S7 [28]. Repeating the protocol 
outlined in Fig. 3a but now as function of 0 with a suf¬ 
ficiently intense pump so that correlations are generated 
with perfect fidelity yields the result shown in Fig. 3g. 
Here the parallel spin alignment from the two modes as 
captured in Fig. 3f can be continuously tuned to an anti¬ 
parallel alignment via (p as shown in Fig. 3h correspond¬ 
ing to an effective phase transition from ferromagnetism 
to anti-ferromagnetism or in other words Jik —Jik 
the Ising model. 

Indeed the pulse sequence shown in Fig. 3a can also 
be executed with the off resonant excitation of the asym¬ 
metric mode in the second step and triggering via the 
symmetric mode as detailed in Fig. S8. Interestingly, the 
data shown in Fig. 3c (and Fig. S8c) indicates that corre¬ 
lations from the two-mode squeezing can be observed via 
the parametric resonances with a pump amplitude that 
is three orders of magnitude smaller than the equivalent 
measurement utilising thermomechanical noise fluctua¬ 
tions as shown in Fig. 2d. This observation suggests an 
alternative approach to identifying mechanical vacuum 
squeezed states via the large signal from their parametric 
resonance instead of their vacuum noise [30| which could 
prove pivotal in detecting a macroscopic all-mechanical 
entangled state [22] . 

The results detailed in Fig. 3 confirm that the requisite 
features of the Ising Hamiltonian can be reproduced in 
the electromechanical domain. However its implementa¬ 
tion with two mechanical modes in this instance is trivial 
and therefore it is instructive to examine the prospects of 
an electromechanical Ising machine composed of a large 
number of spin particles with multiple degrees of cou¬ 
pling. 

Based on the above demonstrations, an electromechan¬ 
ical Ising machine visualised in Fig. 4 would seem feasi¬ 
ble. Here multiple mechanical elements with different fre¬ 
quencies are weakly mechanically coupled to their neigh¬ 
bours. The elements encode spin information via the 
bistable phase of their piezoelectrically activated para¬ 
metric resonance in their fundamental mode, from the 
right clamping point, where this information can also 
be pre-programmed [9]. The key difference here is that 
spin information is stored in each mechanical element 
rather than modes composed of more than one element 
as demonstrated above. On the left clamping point a 
coupling gate electrode is defined which interconnects the 
mechanical spins via the piezoelectrically activated para¬ 
metric down-conversion pump at the sum frequency of 
two elements for instance uji + as shown in Fig. 4. 
Naturally, the reduced mechanical coupling between the 



FIG. 4: The electromechanical Ising machine. A con¬ 
ceptual image of an electromechanical Ising machine based on 
an array of mechanical elements each parametrically resonat¬ 
ing to encode classical spin information (red corresponds to 
spin up and green to spin down) where two-mode squeezing 
is used to create coupling between elements. The parametric 
resonances are piezoelectrically activated and readout via the 
individual gate electrodes (orange) on each mechanical ele¬ 
ment where the element has a natural frequency uJi. The 
global gate electrode (orange) located on the left clamping 
point of all the mechanical elements can enable coupling be¬ 
tween any pair of mechanical spins with the application of a 
sum frequency pump modulation. Utilising FDM, the pump 
can execute multiple degrees of coupling between a large num¬ 
ber of spins potentially enabling the Ising Hamiltonian to be 
solved in a configuration that is beyond the reach of conven¬ 
tional computers. 


elements will require a stronger pump amplitude which 
would be feasible as even weak pumping can efficiently 
couple the mechanical spins as demonstrated in Fig. 3d. 
Each mechanical spin could then be readily coupled to all 
its neighbours by employing a frequency division multi¬ 
plexed (FDM) pump composed of multiple sum frequen¬ 
cies where Fig. 4 explicitly depicts the FDM pump input 
necessary to couple the element to its first six near¬ 
est neighbours. The compact form of this coupling is in 
stark contrast to the optical Ising machine which requires 
delay lines that increase both in number for more spins 
and length for higher order couplings [21]. The FDM 
piezoelectric pump in principle enables a large number 
of mechanical spins to sustain multiple degrees of cou¬ 
pling permitting the electromechanical Ising machine to 
be programmed with problems that are unsolvable with 
conventional computers. Although this platform could 
tackle NP-hard problems 1201 [21], it does not offer speed 
up as it employs classical annealing where thermome¬ 
chanical fluctuations of the mechanical elements drive the 
search in the underlying potential energy landscape for 
the ground state. However tantalisingly if all the mechan¬ 
ical elements are operated in their ground state [TTHTS] 
then quantum tunnelling could be harnessed to increase 
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the speed of search thus making an inexorable case for 
development of such an electromechanical Ising machine. 
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